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Q ' Abstract 

' In this paper we present a new family of discrete sequences having "random hke" uniformly 

^ ' decaying auto-correlation properties. The new class of infinite length sequences are higher order 

00 ■ 

lO ■ chirps constructed using irrational numbers. Exploiting results from the theory of continued 

' fractions and diophantine approximations, we show that the class of sequences so formed has 

' the property that the worst-case auto-correlation cocfhcients for every finite length sequence 

OO , decays at a polynomial rate. These sequences display doppler immunity as well. We also show 

that Toeplitz matrices formed from such sequences satisfy restricted-isometry-property (RIP), a 
, . ■ concept that has played a central role recently in Compressed Sensing applications. Compressed 

' sensing has conventionally dealt with sensing matrices with arbitrary components. Nevertheless, 

such arbitrary sensing matrices are not appropriate for linear system identification and one must 
employ Toeplitz structured sensing matrices. Linear system identification plays a central role in 
a wide variety of applications such as channel estimation for multipath wireless systems as well 
as control system applications. Toeplitz matrices are also desirable on account of their filtering 
structure, which allows for fast implementation together with reduced storage requirements. 
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1 Introduction 



Sequences with special properties are required in a number of applications ranging from communi- 
cation systems, radar systems and system identification. Many of these applications require families 
of sequences with low auto-correlation, cross-correlation and resiliency to doppler offsets. In com- 
munication systems signature sequences with low auto-correlation properties have been employed 
in wireless commmunications (see [14]) to overcome self interference due to multi-path effects and 
co-channel interference due to multi-access communications respectively. 

Tocplitz structured matrices naturally arise in linear system identification, a problem that is 
at the core of many applications ranging from channel estimation in multipath wireless systems to 
model estimation in control applications. In wireless systems the channel coefficients are sparse 
and channel estimation is identical to compressed sensing with Toeplitz structured matrices. In 
parallel in several control applications, signals with low auto-correlation sequences are required 
when a model parametrization is generally unknown and models of increasing complexity are adap- 
tively chosen [16] to characterize the system. This usually requires a persistently exciting input of 
arbitrary order. In addition, an input probing signal that is '"optimal" to any model order is prefer- 
able and this leads us to seek aperiodic input sequences with good non-circular autocorrelation 
properties. It can also arise in the context of system identification in the presence of unmodeled 
dynamics and noise. There is a wide array of literature dealing with system identification in the 
presence of unmodeled dynamics (see for example [13] and references therein). The solution to 
this problem requires designing inputs that can suppress the corruption due to both noise and 
unmodeled errors [20]. While noise can be overcome by applying a persistent input, to suppress 
the effect of unmodeled dynamics requires aperiodic input sequences whose non-circular (or zero- 
padded) autocorrelation uniformly decreases to zero. It turns out [18, 20] that the contribution 
from the unmodeled dynamics in the parametric error is directly related to the rate of decay of the 
worst-case non-circular autocorrelation coefficients. 

In many applications such as multipath wireless systems [3] as well as other mediums such as 
acoustic/RF the system can be modeled as a finite impulse response (FIR) sequence. It turns out 
that in these cases the support of FIR [3] is sparse. Consequently, this leads to the problem of 
reconstruction of sparse FIR systems. 

Sparse reconstruction has been extensively dealt with in the context of Compressed Sensing 
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(CS) [5-7, 9]. CS involves recovering a sparse signal x from linear observations of the form y = Ux. 
Suppose the true signal, x, has fewer than k non-zero components the solution to the £o problem 
recovers this solution if and only if every sub-matrix of U formed by choosing 2k arbitrary columns 
of U has full column rank. Unfortunately, it turns out that the general solution is intractable. In 
[5-7, 9, 19] it is shown that for sufficiently small k, the ii problem is equivalent to solving the io 
problem whenever the sensing matrix U satisfies the so called RIP property. The RIP property 
of sparsity k amounts to the requirement that the singular values of every sub-matrix formed by 
selecting k columns of U is close to one. The question now arises as to what matrices satisfy RIP 
property. Earlier work on CS established RIP property for random constructions and recent work 
has produced deterministic designs based on deterministic sequence constructions [8, 11]. 

Nevertheless, this setup is not directly applicable to the sparse FIR system identification prob- 
lem. Specifically, in this setting the output is the convolution of inputs with the channel coefficients. 
This amounts to taking a product of Toeplitz matrix of the inputs and a vector consisting of FIR 
coefficients as components. Therefore, this calls for RIP property for Toeplitz matrices. Motivated 
by these issues [1] consider random Toeplitz matrices and derive information theoretic limits while 
[3,4] develop RIP properties for random generated Toeplitz matrices. In this paper we design 
deterministic Toeplitz matrices having guaranteed RIP properties. 

To account for compressed sensing in the context of system identification we first present a new 
family of discrete sequences having "random like" uniformly decaying auto-correlation properties. 
The new class of infinite length sequences are higher order chirps constructed using irrational 
numbers. Exploiting results from the theory of continued fractions and diophantine approximations, 
we show that the class of sequences so formed has the property that the worst-case auto-correlation 
coefficients for every finite length sequence decays at a polynomial rate. We also show that Toeplitz 
matrices formed from such sequences satisfy restricted-isometry-property (RIP). Linear system 
identification plays a central role in a wide variety of applications such as channel estimation 
for multipath wireless systems as well as control system applications. Toeplitz matrices are also 
desirable on account of their filtering structure, which allows for fast implementation together with 
reduced storage requirements. Additionally, we show that these sequences are immune to doppler 
offsets. 

The organization of the paper is as follows. In Section 2 we motivate Toeplitz structured 
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matrices in CS and System Identification. We then formulate the problem of CS in terms of 
seeking sequences with autocorrelation properties. This motivates the design of sequences which 
is dealt with in Section 3. The following Section 4 deals with matrix constructions that further 
improve upon the RIP properties. We finally conclude with a brief discussion of practical and 
implementational issues. 



2 Problem Setup 



Linear system identification problems are generally characterized by measured output, y, the input, 
u, Gaussian noise, w, and a linear time-invariant system, G that are related by the following discrete 
time equation. 

y{t) = Gu{t)+wit), t = 0,l,...,n, (1) 

where, the system G is identified by its kernel, {gt}, which is usually assumed to be an element of 
ii, where, ii is identified with the space of bounded analytic functions on the unit disc [20]. The 
objective is to estimate G in some metric norm. Obviously, with finite data arbitrary elements of 
£i cannot be estimated. One typically assumes in this context that the system can be decomposed 
into a finite dimensional model, H and a small residual unmodeled dynamics, A. If the finite 
parametrization happens to be the class of finite impulse response sequences (FIR) of order m we 
obtain in expanded notation, 

FIR Model Residual Error 



ys = Gus + ^5 = ^ 9kUs~k + ^ 9kUs-k +^8 

k=0 k=m+l 

The task is to estimate H = {go, gi, . . . , gm} given that A = {gm+i,gm+2, ■ ■ 
than 7. 

Further expanding into matrix notation we get. 



(2) 

has ii norm smaller 



y2 



Uo 
Ul 





Uo 



91 

92 



+ 



W2 



(3) 



\ Un / I Un-1 Un-2 ■■■ Uq \ \ gn J \ W-n ) 

Now in the absence of noise it is clear that a pulse input of unit-amplitude is sufficient to recover 
the FIR model exactly. On the other hand persistent noise can only be averaged out by a persistent 



input. However, a persistent input will also "excite" the unmodeled error. In particular for the 
above equation, whenever, s > m + 1, the data also contains contributions from the unmodeled 
dynamics. For a periodic input, u, with Us+i = Ug, the data for s = 0, I, 21, . . . can be written as: 



Thus one obtains information only on linear combinations of the g^^s and not the individual coef- 
ficients. It is therefore impossible to determine only the model-coefficients no matter how long the 
input signal and how large the length of the period. Therefore, no matter how large the period, in 
the worst case the unmodeled error will couple with the model-set dynamics. Although we describe 
the simple case of FIR models and residual error, similar decompositions can be obtained for more 
general parameterizations [18,20]. 

We define aperiodic autocorrelations function to characterize the design problem for robust 
identification. 

Definition 1. A n- window aperiodic autocorrelation function (ACF) of a real/complex valued 
sequence, {ut}tez+i is defined as: 



where the asterisk denotes the complex conjugate. 

It turns out that a sufficient condition (see [18] for necessary condition) on input sequences for 
estimation of the optimal finitely parameterized model is that the worst-case aperiodic autocorre- 
lation asymptotically approach zero. 




0<T<n |f^(0)| 

The special case when the unmodeled dynamics, A, is negligible is of importance in many 
applications. Here the impulse response terminates after p coefficients and an input of length n is 
usually applied. Consequently, we get 



s 



ys = ^ 9kUs-k + ws = {go + 9i + ■■■)uo + {gi + gi+i + . . .)ui + . . . + 

k=0 



's 



" s=0 




s 




k=0 



and 



p 



Vs 



^QkUs-k + Ws, s>p 
k=0 
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This corresponds to a "fat" Toeplitz matrix of the form: 



Uo 
Ui 



U = 



Un-l 



Uo 
Ul 



(4) 



Un-l 

where the number of columns in the above matrix is p. 

There is an other situation that arises in the FIR context when the output response is observed 
in steady state. Here, the input ut is assumed to start from time t = —oo and the output is observed 
between time t = p until t = p + n — 1. In this situation we get, 

p 

ys = ^ 9kUs-k + Ws 
k=o 

Organizing the above equation as a matrix for s > p leads to the following Toeplitz matrix: 



U 



Up 
Up+l 



Up-l 



U2 Ul 



Us U2 



U„ 



(5) 



Up+n-l '"n+fe-2 

Note that the structure of the steady state matrix leads to time-dependent auto-correlations, which 
we define below: 

Definition 2. An n- window autocorrelation function (ACF) of a real/complex valued sequence, 
{ut}tez+j is defined as: 

^ t+n 



n 

s=t 

where the asterisk denotes the complex conjugate. 

So far the discussion above has focused on different aspects of system identification. We will 
now relate this to compressed sensing problems. These problems arise when the support of FIR 
coefficients is sparse. As we described in the previous section a central ingredient typically employed 



in establishing that convex programming algorithms lead to sparse recovery is the so called restricted 
isometry property (RIP). We describe the RIP property next. We follow closely the development 
in [2]. 

Let D be the diagonal matrix consisting of ^2 norms of each column of U. Let, Tq be the 
collection of all subsets of {1, 2, . . . , p} of cardinality q. Enumerate the elements of this sub- 
collection as / = 1, 2, . . . , (^) and let 7r(Z) denote the 1th element of Tq. Ut^(^i), D„(i^ denote the 
sub-matrices of U, D respectively obtained by selecting columns with indices in 7r(/). 

Let, Til be the normalized correlation matrix, i.e., 

We next define the RIP property. 

Definition 3. [RIP property] A nxp sensing matrix, U, characterized by the triple {n,p, q) is said 
to have RIP property of order q if there is a number dq G [0, 1) such that every correlation matrix 
satisfies: 

(1 - Sq) < Xnuni^l) < Ame«(S,) < (1 + 5,), V Z = 1, 2, . . . , 

Note that RIP property of order q immediately implies RIP property of smaller order. This 
follows from the observation [2] that for tt (Z) C 7r(Z) we get 

Amin(C^J(;)C/7r(/)) < Ai„in(C^J(,) , Amax ) < Xmax{U^(^i)U^{l)) 

Wc arc now left to establish eigenvalue bounds for the correlation matrices of all orders. First 
note that the coefficients of the correlation matrix arc the autocorrelation coefficients. As remarked 
earlier the autocorrelation coefficients are different for the different matrices. In particular for 
Equation 5 we have, 

; t, s G 7r(Z) G Tq 

t,s 

On the other hand for Equation 4 we get 

; t, s G 7r(/) G Tq 

t,s 

We can now state a sufficient condition for RIP property based on eigenvalue bounds for cor- 
relation matrices of [10]. These bounds are based on a straightforward application of Gersgorin 
theorem. The reader is also referred to Theorem 2 of [3] for a closely related result derived for 
random Toeplitz matrix constructions. 
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C{t,s-t) 



Theorem 1. Suppose U is Toeplitz matrix of Equation 5. Then U satisfies the RIP property of 
order q if 

B= n,« T If ■;-:"< 1 

Now suppose ?7 is a Toephtz matrix as in Equation 4. Then U satisfies the RIP property of order 
gif 

I, ^ Km 

Proof. The proof follows by a straightforward application of Gershgorin's theorem applied to 
Correlation matrices (see Page 60 in [10]). In particular from [10] it follows that the eigenvalues 
of the correlation matrix can be upper and lower bounded by the sum of the off-diagonals. In 
other words, if 

then, 

Amm(Si) > 1 — Ri, Amax(Si) < I + Rt 

Now to ensure every correlation matrix satisfies the RIP property we take the maximum over all 
the elements in Iq. The condition for Equation 4 follows in an identical fashion and is omitted. ■ 

Our first attempt at establishing RIP property will be by quantifying the worst-case autocorre- 
lation coefficient. In particular in Section 3 we will describe sequences whose worst-case autocorre- 
lation coefficients decay at a polynomial rate. In the following Section 4 we will directly deal with 
matrix constructions. 

Definition 4. A bounded sequence, {u{t)}, t = 1, 2, . . . , n + p — 1, is said to have a uniformly 
polynomial decay (PDACF) property if there is a 7 > such that the worst-case autocorrelation 
satisfies: 

Corollary 1. The Toeplitz matrix A of Equation 5 RIP property of order q = 0{n'^) if the sequence 



< Cn~^, yt = l,2, ...,p, 0<t + T <p 



satisfies PDACF property with decay coefficient 7. 

Note that the PDACF property may not hold for arbitrary number of columns n and we make 
this precise later. 




Figure 1: Comparison of Worst-case auto-correlation functions for various lengths of a 2^ — 1 order 
PRBS, sine-sweeps and Higher-order-chirps or fast sine sweeps (FSS). 

In summary a simple means to realize RIP property is to obtain characterizations with low au- 
tocorrelation. It can be argued that for all practical purposes one could possibly take truncations of 
pseudo-random binary sequences (PRBS) sequences with arbitrary long period. However, it turns 
out that truncations of PRBS sequences do not necessarily have small autocorrelations. Unfortu- 
nately the autocorrelation expressions in Theorem 1 involves correlations of truncated sequences. 
In Figure 1 the worst-case auto-correlations of truncated PRBS of period 2^^ — 1 has been plotted. 
As seen such truncated sequences can have poor autocorrelations. Finally, it is perhaps surprising 
that the widely used sinc-sweep, u{t) = exp(iai^), a G M, i = 0, 1, 2, . . ., does not meet the 
requirements either. In [20], we show that, 

limsup max Ir^X"^)! ^ ^/^ 

0<T<n 

The basic reason is that the auto-correlation function turns out to be equal to sin(nTa/2) / (n sin(Ta)) 
and subsequences tj, and rij can be suitably chosen so that the limiting value of Tja mod (27r) 
tends to zero 

We next provide sufficient conditions in terms of worst-case decay of autocorrelation coefficients. 
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Definition 5. A bounded sequence, {u{t)}, t = 1, 2, . . . , n, is said to have an aperiodic uni- 
formly polynomial decay (PDACF) property if there is a 7 > such that the worst-case aperiodic 
autocorrelation satisfies: 

' < Cn"', Vr = 1, 2, n 



Corollary 2. The Toeplitz matrix A of Equation 4 has an RIP property of order q = 0{n'^) if the 
sequence satisfies PDACF property with decay coefficient 

Note that the PDACF property may not hold for arbitrary number of columns n and we make 
this precise later. 



3 Sequence Design 

It turns out that higher-order-chirps (HOC) sequences are PDACF sequences. These are determin- 
istic sequences of the form 

u{t) = exp(z27rat^) 

where a is some irrational number. In applications the real part of the complex sequence will 
actually be applied and to simplify our exposition we consider the complex valued signal here. 

Theorem 2. The third order chirp sequences has the PDACF property for the ACF function as 
defined in Definition 4 with 7 = 0.25 for r < An for all quadratic irrational numbers. 

Remark: Recall that algebraic numbers are those numbers that are roots of any polynomial 
defined over a field of integers. 

Remark: Note that the higher-order chirps are complex valued. However, identical properties 
hold for real-part and imaginary part of these signals. 

The proof is broken down into several steps. In the first step we will upper bound the ACF 
function by means of simpler functions. 

Denote by ||a;|| the distance of x from the closest integer: 

= min \x — A\ (6) 

II II ^g^l I V / 

and 

[x] = argmin^g^lx - A\ 
The following elementary facts follows: 
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Proposition 1. The following hold: 



\\kx\\ = \\k\\x\ 
a;|| 



|sin(27ra;)| G 
Proof. To prove the first assertion we note that: 



2 ' 



[x + k] = [x]+k, y ke Z+ 



Next we proceed as follows: 



\kx — [kx]\ 

\kx - [kx - k[x] + k[x]]\ = \kx - [k{x - [x])] - k[x]\ 
\k{x-[x])-[kix-[x])]\ 



To obtain bounds for the sinusoid function, we observe that: 



lsin(27ra;)| = | sin(27r||a;||)| G [0.5||a;||, 



where the last assertion is a standard result in elementary calculus. 



With these preliminaries we can reduce the expression for ACF in terms of the function || • || as 
in the following lemma. 

Lemma 1. The ACF function, |r"(r)|, r G Z'^, for the HOC satisfies the following inequality for 
any irrational number a G H. 

fe=0 " " 

and for the aperiodic ACF 

\ruir)\^ <-y ^ l<T<n (7) 

k=o 

Proof. First, the ACF function and its aperiodic counterpart can be simplified by straightforward 
algebraic manipulations as follows: 

k=l 
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sin(27ra/cT(n — r)) 



sin(27raA;r) 



C C 1 ^ 
n n j-^ 



sin(27rLA;x) 



sin(27r/ca;) 



, L = n—T, X = ra, 1 < k,T < n 



Similarly, 



fc=l 



sin(27rQ;A;r7T,) 



sin(27raA;r) 



C C 1 V- 
+- = -+^> 

k=l 



s\'n.{2imkx) 



sin(27rfca;) 



, X = Ta, 1 < k,T < p = Xn 



To further simplify these expressions we consider the function, Hl{-): 

3L 



Hl{x) 



l + L\\x\ 



and obtain the following property. 



Proposition 2. The function Hl{x) is a monotonic function over L, i.e., Hj{x) < Hk{x) for j < k. 
Furthermore, it satisfies: 

sin(27rLx) , , 



sin(27ra;) 



Proof. Consider the function, Hl{-): 



Hl{x) 



3L 



1 + L||a;|| 



(8) 



Since, | sin(27rLx)/ sin(27rx)| < L, it follows 

(1 + L\\:r\\) siu(27rLx) 



sin(27ra;) 



<L + L— 



sin(27rx) 



|sin(27rLx)| < 3L 



where for last inequality we have used Proposition 1 



Next, we show that Hl{-) is monotonic in L. Suppose, we have given two integers, n, m, with 
n < m, it follows that: 



m + mn||x|| < n + mn||x| 
The result now follows by inspection. 



m 



< 



n 



1 + m X 1 + n X 



Note that the upper bounds for both cases are identical for A = 1, i.e., p = n. For simplicity 
we consider A = 1 case and comment on what happens when A > 1 later. The main difficulty now 
in establishing the result is that liminfg \\qa\\ = for every real number, i.e., \\ja\\ comes close 
to zero infinitely often. The proof therefore rests on the fact that very few terms in the sequence 

(||ra||, ||2ra||, ||nra|| are close to zero for any phase < r < n. A well known result in 
continued fraction expansion theory [12] provides how closely can rational numbers approximate 
irrational numbers: 
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Theorem 3. If a is a quadratic irrational number(i.e. irrational solution of a quadratic polynomial) 
there exists a constant > such that, 

Furthermore, for every constant e > 0, for almost all irrational numbers a G (0, 1) except on set of 
lebesgue measure zero, the inequality, 

lk«ll < 1 \+e , qe z+ 

glog^+^g 

has only a finite number of solutions. 

Remark: To put things into perspective we point out that Cq, = ^ for the golden ratio, 
a = (1 + V5)/2. 

Unfortunately, this inequality alone is insufficient for deriving the upper bounds for worst-case 
ACF, as shown below: 

1 " 1 1 " 1 

max — > — < — > rr- -/-^ 

o<T<n n ^-^ 1 + n\\kTa\\ n ^ I + 

k=0 k=l ^ kn 

Consequently, we need to rely on more intricate properties of continued fraction theory to develop 
an accurate estimate. 

The above setup suggests that it is convenient to deal with the family of n irrational numbers, 
ra with < r < n. We denote by the symbol /? any irrational number in this family. Specifically, 
we let 

PeSia,n) = {\\Ta)\\; r = l, 2, ...,n} 

Note that by employing Proposition 1 we can further simplify the expression in Equation 7 as 
follows: 

1 " 1 1 " 1 

— > ni n" = — / II , for some (3 G S(a, n) 

n^^l + n\\kTa\\ n ^ 1 + n||/e(||/3||)|| ^ ^ 

Consequently, we need to show that: 

1/2 

<n-^, V/3G<S(a,n) 



(if ' V 



We summarize properties from elementary continued fraction theory for the purpose of completion 
next. These properties and the accompanying notation have been adopted from [17]. 
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Continued Fractions: The continuous fraction expansion of any positive irrational number {3 is 
given by: 

13 = [ao; ai, 02, . . . , aj, . . .] = ao H ^-^ , aj G Z+ 

"2^ — 5~ 

...Qj... 

where for the sake of brevity the first expression is typically used to denote the expansion. In our 
context we have /3 < 1/2 and so oq = 0. By truncating the continuous fraction expansion we obtain 
a sequence of rational numbers called convergents, i.e., 

— ^ = [ao;ai, a2, . . . , ak] 

The numerator and denominator of the convergents satisfy the following recursion: 

Ak+i = ak+iAk + Ak-i, A-i = 1; = ao, k>0 (9) 
Bk+i = ak+iBk + Bk-i, 5_i = 0; So = 1, A; > (10) 

Since, a^s are positive integers the numerator and denominator sequences B^, A^s form a strictly 
monotonically increasing sequence. It follows that: 

Bk = auBk-i + Bk-2 = {atak-i + l)5fe-2 + 5fe-3 > 25fe_2 > 2(^-i)/2 (11) 
where the last equation follows from the fact that > 1. 

Optimal Approximation: The convergents are optimal rational approximants, in that, the kth 
convergent is the best approximation to the irrational number among all rationals which have a 
smaller denominator: 

\\BkP\\ < \\j(3\\, yO<j<Bk, k>l (12) 

For convenience we denote by, 

Dk = BkP - Ak 

The DkS follows a simple a recursion identical to Equations 9, i.e., 

Dk = akDk-i + Dk-2, D_i = -1 (13) 

It follows that the so called even convergents, A2k/B2k, monotonically approach (3 from the right 
while the odd convergents, ^2fc+i/-B2fc+i, approach j3 from the left as shown in Figure 2. It turns 
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H — hH — I 1 h- 



A A A A 

•^2k ■^2k+2 ■^2k+3 •^2k+l 



^2k ^2k+2 ^2k+3 ^2k+l 



Figure 2: Illustration of how even and odd convergents approach an irrational number a 

out that the value of Di^ is bounded from above and below by the size of the denominators given 
by: 

^ < „ < \Dk\ < ^ (14) 



Ostrowski Decomposition: In the so called Ostrowski's representation, each integer, m, is 
expressed as a linear combination of the S^s. It turns out that this representation is unique in the 
following sense. Any integer m G Z'^ with Bp <m < Bp+i can be uniquely decomposed as: 



p 

m 



j=0 

with < Cj+i < Oj+i for j > 1 and < a < ai. Moreover, cj = if Cj+i = flfe+i. This 
decomposition is obtained by first dividing m with the largest denominator. Bp, smaller than m 
and then decomposing the remainder with a corresponding largest denominator and so on. 

Lower Bounds to Rational Approximation: The significance of Ostrowski representation 

follows from upper and lower bounds. First note that due to translation invariance we have, 

p p p p 

w^pw = II XI ^j+i^j^ii = iiX^j+i^j^~X'^j+i^jii = iiX^j+i^jii 

j=0 j=0 j=0 j=0 

The following Lemma from [17] (Lemma 1 Chap. 2) will be useful: 
Lemma 2. 

p 

j=0 

To derive bounds we need to define the type of an integer m. Let, m = '}2^=o'^j+^^j' 
Ostrowski representation for m. The type 7(m) of an integer m is an integer defined as: 

7(m) = min{j | cj+i ^ 0} 
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We also define 7* to be the maximum possible type in Z„ = {1, 2, . . . , n}, i.e., 

7* = max{7(m) \ m = 1, . . . , n} 
Proposition 3. The maximum type, 7*, in the set Z„ is smaller than 21og2(2?7.). 
Proof. If m e Zji is of type 7(m), then, from the ostrowski decomposition it follows that 



7(m)-l 

m > B^(^) > 2 2 



2 log2 m + 1 > 7(m) 



where, wc have used Equation 11 for the second inequality. Thus, the largest possible type in Zn 
is A„ = 2 log2 n. Since, type is a positive integer ranging from to 2 log2 n, the result follows by 
inspection. ■ 



The type of the integer, q, signifies how small \\qP\\ can get. We have the following result from 
[17] (Theorem 1 Chap. 2) that precisely characterizes this connection: 

Theorem 4. Let < /3 < 1/2 and the integer, m, have type j{m) > 1. Then, 



\m, 



i=o 



The above results lead to the following precise characterization of how type is associated with 
\\mP\\: 

Lemma 3. Let < P < 1/2 be an irrational number and m > 1 a positive integer. Then, 



mB > < 



min ((ci - l)|L>o| + (02 - C2)\Di\, (oi - ci)|L>o| + C2I-D1I) 
> min 



2B2 ' 2Bi 2B2 J 



2-B7(m)+l 



+ 



2-B7(m)+2 



7(m) = 



7; - 1 > 7(m) > 1 



7("^) = 7n > 1 
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Proof. The case of j{m) = 7* is a direct consequence of Theorem 4. To estabhsh the second case 

we proceed by noting that for 7(m) > 1 we have \\m(3\\ = \ Yl^=o^j+^^j\ from the fact that 

the even and odd D'/^s have opposite sign we get: 
p 

I ^ = \Cj(^m)+lD-y{m) + '^■y{m)+2D^(^rn)+l + C-y(m)+3-^7(m)+2 

j=0 

+ • • • + C^{m)+2jD^(^m)+2j-l + C^(m)+2j+lD.yl^rn)+2j + • • • | 
^ |c7(m)+l-C^7(m) + C^{m)+2D^(^rn)+l + 'J7(m)+4-C'7(m)+3 + '^7(m)+6-C'7(m)+5 + ■ ■ ■ | 
^ \Cj{m)+lDj(rn) + a7{m)+2-D^(m)+l + (c7(m)+2 — ('■j(^m)+2)Dj(^m)+l + a7(m)+4-D^(m)+3 

+«'7(m)+6-f^7(m)+5 + • • • I 

The second inequaUty follows from the fact that we have removed all the terms with the same sign 
as terms of the form -D^(r„)+2j) while taking the maximum possible values for terms 

of opposite sign. The third inequality is just a restatement of the third term. Next we utilize 

Equation 13 to make the substitution that — Dk^2 = CLkDk-i and get, 

p 

I ^ Cj+i-DjI > |c^(m)+l-^7(m) + (c7(m)+2 " 0,-y(m)+2)Dj(^rn)+l + «'7(m)+2-f^7(m)+l 
j=0 

+a7(m)+4-C^7(m)+3 + 'l7(m)+6-C^7(m)+5 + • • • | 
^ |c7(m)+l-^7(m) + (c7(m)+2 ~ 0'-y{m)+2)Dj(^m)+l 

+ {Dj(m)+2 - -D7(m)) + iDj(m)+4 " D^{m)+2) + • • ■ | 
^ KSM+I ~ l)-^7(m) + (c7(m)+2 ~ 'i7(m)+2)-^7(m)+l I 

Finally, we know from Equation 14 that 

1 1 

\D^(rn) \ > — > 



Next observing that -D-y(m) and are of opposite signs the result follows. 

Finally, to establish the first case, i.e., 7(m) = we note that a lower bound for | ^^^=0 
follows identically as the case for 7(m) > 1. We next derive an upper bound by excluding all terms 
of opposite sign (except for the term containing D-\) as Dq and include the maximum possible terms 
of the same sign. Through algebraic manipulations as for the lower bound and noting that the even 
D'^s are of the same sign we get: 

I ^ Cj+iDj| < \ciDq + C2Dx + a3D2 + 0^04 + . . . | < \ciDo + (c2 - l)Di\ 
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Next, substituting Equation 13 for Di = uiDq — 1 we get: 

In 

\^Cj+iDj\ < |1 + (ci - ai)Do + C2D1I 

i=i 

Now it follows that the latter term is strictly less than 1 by noting that Dq = 13 > and Di is 
negative and: 

|1 + (ci - ai)Do + C2Di\ < |1 + ((ai - 1) - ai)Do\ < |1 - Do| < 1 

We also note that D2 = 02 -Di + Dq > 0, and so, 1 + (ci — ai)Do + C2-D1 = ciDq + (c2 — l)Di > 
ciDq + (02 — l)Di > (ci — 1)Dq — Di + D2 > 0. Thus we have a upper and lower bound for 
I Yl]=i Cj+i-DjI that are both positive and smaller than 1, i.e., 

< (ci - l)Do + (C2 - a2)Di < I ^ Cj+iDj\ < 1 + (ci - ai)£>o + C2-D1 < 1 

Now noting that Dq is positive and Di is negative establishes the first inequality for 7(m) = 0. Also 
substituting the lower bounds of Equation 14 it follows that the second inequality for j{m) = is 
also satisfied. ■ 

To complete the above set of results we need a bound for the situation when 7(m) = 7* = 0. 
We have the following corollary: 

Corollary 3. Let < /? < 1/2 be an irrational number and m > 1 a positive integer and 7(m) = 
7;^ = 0. Then, 

\\mP\\ > min ((ci - l)|Do|, (ai - ci)|Do|) 

We are now ready to prove the main theorem by a combination of the above well-known prop- 
erties in continued fraction theory. The main outline of the proof is as follows. Lemma 3 and 
Corollary 3 points to the fact that the type of a number controls the value of || ma ||. This motivates 
partitioning of the set Z„ = {1, 2, . . . , j, . . . , n} based on its type and computing the contributions 
for each type. This leads us to defining the following sets: 

•^1,0 = {m G Zn I 7(Tn) = I, c^(rn)+i = c, in the Ostrowski expansion for m} 
•^1,0,(1 = {m € Zn I j{m) = I, c^(^)_|_i = c, c^(^rn)+2 = d in the Ostrowski expansion for m} 

We show in the sequel that the cardinality decreases exponentially with the type. 
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Lemma 4. The cardinality of Ai^c in the set Z„ is given by: 



#Ai,c < { 



2 



Kin 



1 l = % 



n 



and the cardinahty of the set, Ai^ci,c2 is bounded by the cardinahty of the set, Ai+i,c2- 
Proof. The proof requires the following proposition. 

Proposition 4. Given, two positive integers, q, r and their corresponding Ostrowski representa- 
tions, {cj{q)}, {cj{r)}, it follows that: 

q <r <J=^ 31 e Z'^ such that < Q(r); Cj{q) < Cj{r), Vj > / 

Proof. (<^= case) Suppose there is an I satisfying the hypothesis. It then follows that, 

p p 
r-q = ^Cj+i{r)Bj-^Cj+i{q)Bj 

j=0 j=0 

1-2 p 

= - Cj+i{q))Bj + ^ {cj+i{r) - Cj+i{q))Bj 

3=0 3=l-i 
1-2 

> -^Cj+i{q)Bj+Bi_i 

3=0 

It remains to show that the last term is positive. We do this by induction. Clearly, the hypothesis 
is true for / = 1. Suppose, the induction hypothesis is true for / = k, then for / = A; + 1, we have: 



(ajk+i - l)Bk + Y!j=l Cj+iBj, Cfe 7^ 
o-k+iBk + Z]j=i Cj+i-Bj, Cfc = 

where, the inequalities follow from induction hypothesis and Equation 9. The RHS corresponds to 
all the admissible CjS and the proof follows. 

case) This follows by contradiction and reversing the previous arguments. ■ 

Proof. (Lemma 4) Any integer q G Ak,c has an Ostrowski decomposition: 

q = ^Cj+i{q)Bj = (0, . . . , 0,Ck+i{q), Ck+2{q), ■ ■ ■,c^*+i{q)), Cfe+1 = c 

j=k 
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where Ck+is are arbitrary integers constrained only by property (C). Consequently it follows from 
Proposition 4, for any other integer, p G Ak c, P > Q that, 



31 > k + 2 such that q(p) > ci{q), Sz Cj{p) > Cj{q), V j > I 



This implies that. 



i-i 



p-q = ^ {cj{p) - Cj{q))Bj_i = ^ {cj{p) - Cj{q))Bj_i> Bi_i - ^ Cj{q)Bj_i{lQ) 
j=k+l j=k+2 j=k+2 

l-l 1-2 

> ai^iBi_2- ^ Cj{q)Bj^i> Bi_2- ^ Cj{q)Bj^i 

j=k+2 j=k+2 

> . ..> -Bfc+2 — Cfe+2-Bfe+i > {ak+2 — Ck+2)Bk+i > -Bfc+i 

This means that there can only be one term belonging to Ak,c for any sequential set of Bk+i 
integers. Now n can be written as: 



n 



n 



B 



k+l 



Bk+i + r, < r < Bk+i 



Therefore, the remainder, r, terms can contain atmost one term. This implies. 



n 



Bi 



+ 1 



Now, since 7* is largest possible type in Z„ we see that -B-y*+i must be larger than n. This implies 
that. 

For all other types, < j{m) < 7* — 1 we have B^{m) < n, which implies, 



#Alc < 



n 



Bi 



l+i 



+ 1 < 2 



n 



, < 7(m) < 7; - 1 



Finally, in order to compute the bounds for i^Ak,c.d we observe that the first two Ostrowski expan- 
sion coefficients for any two integers p, q & Ak,c,d are identical. This implies. 



p-q= {cjip) - Cj{q))Bj_i = ^ (cjip) - Cj{q))Bj_i 

j=k+l j=k+3 
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and the rest of proof follows as in Equation 16. ■ 
Proof. (Theorem 1) The proof of Theorem 1 follows by combining lemmas 4 3. 

l + n\\kB\\ ~ n ^ ^ 1 + n\\m3\\ 

_ 1 ^ 1 2: V V 1 

n ^ , l + n||m/?|| n ^ ^ l + n||m/?|| 

1^1 

Next, we simplify each of the three terms. First we compute the contribution for the second term: 

\ — ^ \ — ^ \ — ^ \ — ^ \ — ^ \ — ^ 

^ l + n||m/?|| ~ ^ \ ^- n\\mB\\^ n ^ ^ l + n\\mB\\ 

(a) r^i';^iCard(A,c) , r^'""^"'Card(A,i,c) 



Card {Ak,c) 1 v;^ v;^ 
f-f ^ 1 + n f-f ^ T 



fc+2 



fc=l c=2 -^^'"-Bfe+i fc=l c=l "T 

(d) i^J^i 1^4^^ 47* log(n) (e) 4(log(n 

n ^-^ n n n 

fe=i fe=i 

The inequality (a) follows from bounds in Lemma 3. We utilize the fact that 1/(1 +x+y) < l/(l+x) 
for positive x and y to split the sum in two parts with terms belonging to Ak,c and A,i,c sets. 
Inequality (b) follows from Lemma 4. Inequality (e) follows from upper bound for 7* derived in 
Proposition 3. 

Now to simplify the first term of Equation 17 we note that: 

11 11 

< 



l + n\\mP\\ i + + 1 + + i + n%^ + n^ 

(18) 

Bounds for the two terms on the right now follows in an identical fashion as the steps for the first 
term derived above. 
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For the last term in Equation 17 we have two cases to consider: j{m) = 7^ > 1 and j{m) = 
7* = 0. For both cases from Lemma 4 the terms of type 7* are of the form Bj*, 2Bj*, ... , d^^B-^* of 
which B^* is the first term in Zn- Now if, a, is a quadratic algebraic number and /? = ra G 5(n, a), 
we know that 



ID II \ 

\B^*Ta\\ > — 



(19) 



Since, the last term, cB^* of type 7* has to be smaller than n, we have, 



n 
B^* 



Since, P = ra for some < r < n, we have, 

F ^ r 

' 1 + n||mra| 



< 



-E — ^^-E^T- = -E 



1 / B * 
< -min(4,--^log(l + d*)) 



< — min 
n 



n 



B^ 



log 1 + 



in 



n 



In 



< 



log^ 



n 



(20) 



Now for the case when j{m) = 7* = we invoke Corollary 3 and follow along the lines of the 
proof for 7(m) = of Equation 18 by first splitting it into two terms, i.e., 
1 1 



1 + n||m/3|| 



< 



1 + nmin ((ci - l)\Do\ + (03 - C2)|-Di|, (ai - ci)|L>o| + C2|£'i|) 

1 1 

+ 



l + n(ai -ci)|L>o| 1 + n(ci - l)|L>o| 
We now follow along the lines of Equation 20 by first substituting for Dq given by Equation 19 
and following along the lines of Equation 20. To complete the proof we need to ensure that the 
PDACF property continues to hold for < fc, r < p = An, with A > 1. To do this we again analyze 
Equation 17. The first two terms depends only on the cardinality of the set Ak,c- Following along 
the lines of Lemma 4 we see that. 



2A 



-Bi+i 
1 



Kin 
l = 7n 



For the third term we observe that. 



E 

'An 



log°-^(n 
(ACa)0-5 



1 + ?7.||mTa| 



< 



n 
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4 Matrices with RIP Property 

In this section we quantify the RIP property for matrices described in Section 2. Consider the 

Toeplitz constructions in Equation 5. It is clear from the sequence designs in the previous section 
togctlicr with Theorem 2 that the RIP property (sec Definition 3) of order (An, n, n IS 
satisfied. However, this bound is quite loose and we establish substantial improvement over this 
bound. We have the following result: 

Theorem 5. Consider the Toeplitz construction of Equation 5 with the elements generated by the 
HOC with a equal to the golden ratio. It follows that for sufficiently large n, the RIP property of 
order (An, n, n^/^/(-\/Alog(n))) is satisfied for this matrix construction. 

The idea is that to establish RIP property of order /c, by Theorem 1 we only need to show that 
for all subsets, /fc C {1, 2, . . . , An} of size fc, 

E 1^(^)1 < E f- E — I — < 1 (21) 

T&Ik T&Ik \ m=l / 

Consequently an autocorrelation decay of order n^'^/^ is only required on "average". Now there 
are contributions of three terms as seen from Equation 17 and the first two terms already decay at 
this rate. Therefore, we are left to analyze the third term, namely, terms that belong to the largest 
type, 7* for each r. These are the only terms that contribute towards the slow decay. Therefore, 
we are left to establish that the number of terms with slow decay are relatively small, which we 
present next. 

Next we decompose the RHS into three terms as in Equation 17. We note that the first two 
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terms are oi 0{n ^). Therefore we can write: 



„2 I 1 + n||i 



1/2 



< 



rSXfc \m=l 



mra 



1 ^ [ A A 1 



1/2 



mra I 



< V^max 



( - V 1 ^ 

V n' n l + n||mTa|| I 



l/2\ 



< V2 



max 



V 



' 1 ^ 1 



V 



^ l + nj||S^*ra| 



First note that by hypothesis we have, \Ik \ < /ri^r ^^'^ hence the first expression in the maximum 
is smaller than one. We are now left to compute the contribution from terms of the largest type. 
Let B^* be the first term of the largest type. Note that B^* depends on r but we suppress this 
dependence to simplify notation. Let 



<^w= E 1^ -iio — ^,^ = 1, 2, 

1 + ni\\B^*Ta\\ 



, n 



We have the following lemma that provides decay bounds for ^(r): 
Lemma 5. Suppose r is such that Bj*t > n^^ and ||B-y.TQ;|| = n"^ for some rji, r/2 G [1, 2] then, 

<^(r) < n'"'"^{i-''i'^2-2} iQg(^) (22) 

for 1 < r]2 < r]i and zero if r]2 > rji. Furthermore, if r]2 < ^ then, 

A 



CnU 



Proof. The last part is immediate and follows from: 



I I 

m 



1 ^ 1 

n ^ 1 + n7||B-y*rQ;| 

7=1 ' -^11 7n I 



< 



Alog(ra) 
Can 



To establish the rest of the statements we note the fact that a is a quadratic irrational. This 
imposes the fact that if B^*t = n^^ then ||iJ-y.ra|| > Can~^^ . Next note that L;^tJ — and 
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n II 5^* rail > nP^ ^. Consequently, when r/i — 1 > 2 — r/2 we use the fact that 



< 1 



1 + jn\\B^*Ta\\ 
On the other hand when 772 — 1 < 2 — r/i we have 

1 1 

< 



1 + jn||i?^*Ta|| jn^ '^^ 

■ 

The above lemma immediately implies that if r is such that r/2 < 1 and Ik < '^fx'^ autocor- 
relation for this shift terms, r, will be negligible. 

Note that the minimum in the RHS of Equation 22 is achieved for r/i = 772 = 1.5. For this 

value wc get 0(t) < \og{n) / y/n. To understand the contributions for different r we partition it into 
different decay rate regions. To this end let, 0.5 = < ^2 < • • • < = 1 be a partition of the 
interval [0.5, 1] with 5 = 6j+i — 6j and L a large fixed integer independent of n. Let 

5,- = {r I (/)(r) G [n~'\ n-^^-^]}, r = 1, 2, . . . , n 

The cardinality of these sets is given in the following lemma: 

Lemma 6. Let a be the golden ratio, i.e., a = (1 + \/5)/2 then for sufficiently large n. 
Proof. Let 

(pir) G [n~^^+i, n"^^] (23) 

It follows from Lemma 5 that. 

Noting that n''^ = -67* t and n^''^'^ = ||i3^* ra || and the fact that a is a quadratic irrational implies 
that ||5^*ra|| > Can~'^^. By substituting for rji and 772 we get: 

||B^.ra|| G [7i-(i+''^+i\7i-(2-''^)] (24) 

In addition we have. 

In ' — 
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Motivated by Equation 24 we denote by 

Cj = {me Z+ \ \\ma\\ G n'^^"^^)], m < n^+^^+i} 

By construction note that for r G Sj and associated largest type B^* the product B^*t G Cj. 
Thus our proof now relies on estimating the cardinality of Cj and then computing the number of 
divisors r for each element of Cj. In other words, 

< \Cj \ max{ # divisors for m e Cj} (25) 

Let Fk, A; = 0, 1, ... be the convergents of the golden ratio. It is well known that the convergents 
satisfy a linear second order recursion 

Fk+i =Fk + Fk-i, Fo = Fi = l (26) 

The continued fraction expansion of the golden ratio is: 

a= [1;1, 1, 1, ...] 

A closed form expression for Fk is given by 

F. = (27) 

where a = (1 + \/5)/2. Prom Lemma 3 it follows that if m G Cj then the type of m in the Ostrowski 
representation (i.e. the smallest element) must satisfy: 

F,im)&[[n^'^'^+'\n('-'^^]] (28) 

Now the largest element in Cj is clearly m = ri^^^i . Consequently, the number of different conver- 
gents, s in the Ostrowski representation of the set Cj can be obtained by employing Equation 27 
and Equation 28: 

(2g,+<5-l)log(n) ^ ^ 
log (a) 

A crude upper bound for the cardinality is clearly 2* since the set of all expansions can be associated 
with s binary digit expansions. However, from Equations 26, 15 we see that no two consecutive 
convergents can be present in any expansion. Therefore, the cardinality of the set can be refined 
and given by the following expression: 

[_s/2\ 



<y ('~^^\< max 2LV2J(-+i^2H) 
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where i?2(-) is the binary entropy function. The RHS is obtained through StirUng's approximation. 
Further simpUfication yields: 

/ ) ^ 2°-^^* = 2°-^5(2^J+''~^)'°s(")/log(a) ^ ^0.78(26»j-|-5-l) 

Finally, to compute the cardinality of Sj we see from Equation 25 we need to determine the number 
of divisors for each m e Cj. This is a classical result going back to Ramanujan [15]: 

Lemma 7. For an arbitrary e > 0, there exists a positive integer mo such that the number of 
divisors is smaller than 2(i+^)'°gM/iogiog(m) f^j. ^ > j^^^. 

Consequently, the cardinality of 5^ for sufficiently large n is given by: 

\Sj\ < 2(l+^)'°gW/loglog(n)^0.78(2e,+5-l) < ^o. 78(26 j+5~l)+v ^29) 

where v is an arbitrary small number for sufficiently large n. I 

Now to complete the proof of the main theorem we determine the set of partitions for which 
the contributions for r G {1, 2, . . . , n} approaches zero, i.e. 

min^ H - mi^X] \^3\^~^' — ' ^ 

j>k TESj j>k 

Upon direct substitution it follows that this is ensured whenever, 0.78(2^j — 1) < Oj/2 (ignoring 
the small 6 and v terms). Upon evaluation we infer that this holds for all 9j < 3/4. We now let 
imin be the minimum j such that 6j < 3/4. Consequently, we are left with 

Jmin 

TGXfc J = lT€lknSj 3>3ininTeSj 

where e > can be chosen to be an arbitrary small positive number for sufficiently large n. Note 
that Oj < 3/8 implies ^(r) < n'^/^. Therefore, if 

Ik < T—r^ 
log(n) 

the RHS can be made arbitrarily small as well and the result follows. I 
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We observe that the golden ratio, a, provides a moderate improvement in RIP property. We 
point to Equation 29 as one of the principle reasons. This Equation shows that number of large 
type of order m = n'^ betweeen two (large) numbers mi = n'^^ and m2 = n'^^ scales polynomially 
with n. This can be attributed to the logarithmic scaling (log{n)) of the number of convergents 
between mi and 1712- Consequently, a question that arises is whether one could control this scal- 
ing by choosing a different quadratic irrational number. One possibility is to choose recursions 
(see Equation 9) that lead to fewer convergents. This can be accomplished by choosing periodic 
continued fraction expansions of the form: a = [ao;ai] where ai is large. However, fewer conver- 
gents does not result in decreasing the number of integers of large type since now we can admit 
B*^Ta, 2B*^Ta, ... , aiB*^Ta as alternatives. Therefore, we conjecture that this is nearly the best 
result one could hope for with quadratic irrational numbers. 

5 Discussion 

We first point out that the RIP order derived in the previous section appears to be small relative 
to what can be obtained with unstructured constructions [11]. We believe that this looseness is 
inherently a consequence of the bounding technique utilized here. To further confirm our belief we 
performed a Monte Carlo simulation to determine the condition number (ratio of the maximum 
to the minimum singular value) for various choices of sparsity levels for Toeplitz matrix of Equa- 
tion 5. In particular we chose the {n,p,q) parameters as follows: We let n denote the number of 
measurements; we set the number of variables p = 2n and the sparsity level q = n/5. The results 
are described in Figure 3. 

We next briefly discuss various implcmcntational and practical issues here. 

First, note that we have described properties for complex-valued signals. Nevertheless it is 
straightforward to check that identical properties hold for real and imaginary parts of the higher- 
order-chirp(HOC) signals. A second claim that can be immediately verified is one of doppler 
resilience. This is because the autocorrelation properties of HOCs are not affected by constant 
frequency shifts. A third aspect is that our construction of Toeplitz matrices could be employed in 
sequential processing. This is because each time a new row is added in Equation 5 or Equation 4 
the new matrix formed by concatenation of the previous matrix with the new row still preserves the 
RIP property. Finally these sequences and matrices can be generated with relatively little memory. 
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100 150 
Number of Measurements 



250 



Figure 3: The condition number (ratio of maximum to minimum singular values) is plotted as 
a function of number of measurements. 5000 Monte Carlo simulations were performed with the 
number of measurements ranging from 50 to 500 for Toeplitz construction of Equation 5. The 
number of variables p was scaled as twice the number of measurements. The sparsity level was 
scaled as 20% of the number of measurements. 
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To see this note that since the autocorrelations decay at a polynomial rate, n~'^, for a toeplitz 
matrix of n rows, we only need to ensure that our approximation does not have larger error than 
this magnitude. Again continued fraction expansion can be employed for this task. For instance, 
consider the HOC 

ut = exp(— j27rat^) = exp(— j27r||ai^ ||) 

For the golden ratio, a = (1 + \/5)/2 we can compute the above expression to any order of 
approximation conveniently. This is because the convergents of a are given by the Fibonacci series, 
Ffc, which can be computed through simple recursion (see Equation 26). Furthermore, to compute 

a rational approximation to qa for any q ^ Z we can proceed as follows. Let Fq be any term in the 
Fibonacci series for which g is a factor (note that there are infinitely many such terms for every 
whole number q [21]). Then 



mm a — — < — = — =^ mm \qa — —\< — =^ qa G 



OP _^ 



V5FI 



q \orq f J. q \ or q 

Since Fq can be chosen to be sufficiently large we can get arbitrarily good approximations through 
employing Fibonacci series. 
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